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We study the interplay of the ferromagnetic (FM) state and the p-wave superconducting (SC) 
state observed in several materials such as UCoGe and URhGe in a totally nonperturbative manner. 
To this end, we introduce a lattice Ginzburg-Landau model that is a genuine generalization of 
the phenomenological Ginzburg-Landau theory proposed previously in the continuum and also a 
counterpart of the lattice gauge-Higgs model for the s-wave SC transition, and study it numerically 
by Monte-Carlo simulations. The obtained phase diagram has qualitatively the same structure 
as that of UCoGe in the region where the two transition temperatrures satisfy Tpm > Tsc- For 
Tfm/Tsc < 0.7, we find that the coexisting region of FM and SC orders appears only near the 
surface of the lattice, which describes an inhomogeneous FMSC coexisting state. 



I. INTRODUCTION 

In the last decade, superconducting (SC) materials 
coexisting with ferromagnetic (FM) long-range orders 
have been found and of intensive interest. In UGe^ 
and URhGe 2 -, a SC state appears only within the FM 
state in the pressure-temperature (P-T) phase diagram, 
whereas in UCoGe^, the SC state exists both in the FM 
and paramagnetic states^. Phenomenological models of 
FMSC materials were proposed 5 ^ soon after their dis- 
covery. The most important observation in those studies 
is that the FMSC state is a spin-triplet p-wave state of 
electron pairsi. 

In the present paper, we propose a lattice model for 
describing FMSC materials and investigate it by nu- 
merical methods^. The model contains a vector poten- 
tial (i.e., gauge field) to describe a FM order parame- 
ter (magnetization) and also a SC order parameter, i.e., 
Cooper-pair field for the p-wave SC state. These two 
physical variables couple with each other as the Cooper 
pair bears the electric charge 2e. Finite magnetization in- 
side the material tends to induce vortices of the Cooper- 
pair field and destabilize SC state. In this sense, the 
FMSC state is a result of frustration. 

Introduction of the spatial lattice in the present model 
has several advantages over the Ginzburg-Landau (GL) 
theory in the continuum space 5 -^; (i) it reflects the lat- 
tice structure of the real materials, (ii) it works as a rc- 
guralization of vortex configurations because the energy 
of these topological excitations becomes finite, and (hi) 
it allows us to make a nonperturbative study by means 
of Monte-Carlo simulations in which contributions from 
all the field configurations including topologically non- 
trivial excitations are taken into account, and so the 
obtained results are quite reliable. In this sense, the 
present study is complementary to the previous analyt- 
ical studies employing perturbative and mean-field like 
approximation^ 5 '" ' ' 



The present paper is organized as follows. In Sec. II, we 
introduce the lattice GL model that is derived the pre- 
viously propose GL theory in the continuum. Detailed 
discussion on the physical properties of the model is also 
given there. In Sect. II. A, we present a brief review on the 
lattice gauge model for the SC state. This review may be 
useful to make the present paper readable for condensed 
matter physicists who are not familiar with the models of 
SC state on the lattice. In Sec. Ill, we present results of 
the numerical study. The phase structure of the model 
is clarified by calculating specific heat, FM correlation 
function, shielding mass of magnetic field, etc. We also 
study behavior of vortices in a constant magnetic field in 
the present model in order to obtain an intuitive picture 
of the Meissner state. Section IV is devoted for conclu- 



II. LATTICE MODEL FOR FMSC STATE 
A. Lattice gauge-Higgs model for SC transition 

In this subsection, we review a typical lattice model 
to describe the conventional SC transition, which is 
called the U(l) gauge-Higgs model (or the Abelian Higgs 
model). Reader who is familiar with this subject can skip 
this subsection and go to Sec. II. B. 

Let us start with the GL theory of s-wave SC state 
in the three-dimensional (3D) continuum space. Its free- 
energy density is given by 
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where ip is the complex scalar field for s-wave Cooper 
pairs, A em is the vector potential (x2e) for fluctuat- 
ing electro-magnetic field, is the covariant derivative 
in the /Lt-th direction (fi = 1,2,3), T® is the bare criti- 
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cal temperature(T) of SC transition, e is the elementary 
charge, a and A are positive T-independent parameters. 

The lattice-field model corresponding to the GL theory 
(|2.1[) is defined by giving its free energy (or the action 
including the inverse temperature) F as follows; 

F = ~ Yl (VS-B^Vfe + cc) +F,, + F A , (2.2) 

X.fl 

F^ =^(a|^| 2 + A|^| 4 ), (2.3) 

X 

U xlt = exp(^), 

where x is the site of the 3D lattice, ip x is a complex 
SC order-parameter field defined on the site x, and we 
use fi also as the unit vector in the /x-th direction. ij) x 
is sometimes called Higgs field because it is a complex 
scalar field. A™ 1 is a real electro-magnetic field put on 
the link (x,x + fj). Fa is the free energy of the electro- 
magnetic field and has the following two versions. One is 
the compact version, 

Fa = - Y ( u w u £+v,n u s+w u *n + c - c -) ' 

A™e(-w,n), (2.4) 
The other is the noncompact version, 

A%e (-00,00), (2.5) 

where V M is the lattice difference operator such that 

Vp/a; = fx+n ~ fx- (2-6) 

These two Fa are distinguished by having the periodicity 
under A xfl -> A xll + 2n or not. We note that F of (|2.2|) 
is invariant under the local U(l) gauge transformation, 

ip x -> = exp(iA a: )V'x, 

40m . /cm _ Mm 1 \ _ \ 

n xy, ^ ■"■Xfj. ~ -"-XH 1 ^x+p, ^x, 

U Xf i -t U' Xfi = exp(iA x+M )f7 :CjU exp(-iA x ) ) (2.7) 

where A^ is a site-dependent real variable. 

Let us first consider the pure gauge system described 
by the energy F A alone. In the case in which fluctua- 
tions of the vector potential A™ are small, the above 
two versions belong to the same universality class, i.e., 
the system is in the Coulomb phase. This is because Fa 
of (|2.4j) approaches to Eq. (|2.5p due to the relation, 

U xv U X j rl/ .uJ^x-^-fiji/Uxfj, ~\~ C.C. = IcOsQxfxv, 
ft — X7 A cm — X7 4 cm 

COS0 X/1) , ~ 1 - \^ x ^v for smau Qxn»- ( 2 - 8 ) 

For large fluctuations of vector potentials, the compact 
version generally allows topologically nontrivial excita- 
tions like magnetic monopoles and may exhibit another 



phase called the confinement phase, which is not allowed 
in the noncompact version. 

Next we consider the case in which the magnetic field 
is switched off by setting A™ 1 = 0. Then the system (|2.2j) 
is reduced to the |(/)| 4 theory. In the 3D |(/)| 4 theory, there 
exists a second-order phase transition accompanying the 
spontaneous symmetry breaking of the global U(l) sym- 
metry under the phase rotation, ip x — > exp(i6)i() x . This 
broken phase is called Higgs phase and corresponds to the 
SC phase. In the limit of A — > 00 with the ratio er/A kept 
fixed to a finite negative value, the system reduces to the 
so called XY model defined with \ip x \ = ^—<j/2\, which 
is well known to exhibit a second-order transition as K 
is varied. This limit corresponds to the London limit of 
the SC (or superfluidity), and the phase transition in the 
\4>\ A theory for large A and that of the XY model belong 
to the same universality class. Even in this simplified 
model, the detailed critical behavior at the phase tran- 
sition is different from that described by the mean-field 
theory (MFT). 

Finally, let us turn on the vector potential A™. In the 
continuum, it is shown^ that the second-order transition 
for ^4 om = is changed to a first-order one as the GL 
parameter k oc A/e 2 is decreased. On the lattice, the 
compact version of the system is studied in RefJ^ and it is 
verified that the phase transition takes place as one varies 
K with fixed k. More precisely, the phase transition is 
of first order for small k and becomes second order for 
large n as in the model defined in the continuum. Here 
we should note that the study in Refi 1 ^ shows that the 
critical behavior of the second-order transition near the 
London limit is not in the same universality class as the 
3D XY model due to the gauge-field fluctuations. 

The noncompact lattice version is also studied in the 
London limit in Ref^. The corresponding energy is ob- 
tained from Eq. (l2.2p by making the replacement, 

Y Y {Vx+^Ux^x + C.C.) 

-+ y E (e-^'+'U^e** + cc) , (2.9) 

where tp x is the phase of the Cooper-pair (Higgs) field 
ip x , and neglecting F$ in Eq. (|2.3l) because it becomes a 
constant. In Refi- 3 this U(l) gauge-Higgs model on the 
four dimensional lattice is studied for large 1/e 2 , which 
exhibits a second-order transition as K is varied. This is 
consistent with the fact that, in the limit of 1/e 2 — > 00, 
the system reduces to the four-dimensional XY model. 
This phase transition is induced by the condensation of 
vortex excitations in ip x . In the terminology of XY spin 
model, the XY spin S x = (cos (p x , sin ip x Y has a definite 
amplitude |<Sa.| = 1, whereas the disorder phase (S x ) = 
is possible due to the strong fluctuations of its angle <p x ■ 
Then it is useful to rewrite the system in terms of topo- 
logical excitations such as vortices and monopoles. In the 
London limit of Eq. (12. 2|) , the duality transformation can 
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be perfome d 13 i , and the free energy is expressed by the 
integer- valued vortex line-element variables J 2/i and the 
integer- valued monopole-density varaiables Q s = V M J Sfl 
sitting on the dual lattice (x denotes its site) as 

x,y \ / 

„ exp(— mnr) Q 

D £jS ~ — ^ ^, r = x - y , ml = 2Ke 2 , 2.10 

m^r 

where D Sj y is the 3D lattice Green's function with mass 
to. As explained in introduction, there appear no singu- 
larities in F v (|2~TU|) . 

For the noncompact version of the gauge system, it is 
shown in RefJ^ that no monopoles exist Q s = and only 
the closed vortex loops that satisfy Ylu^iJ-^f ~ u are 
allowed as expected. In the Coulomb phase with lower 
K these vortex loops condense while in the Higgs phase 
with higher K vortex loops are suppressed. 

On the other hand, for the compact version^, open 
vortex strings may appear and a monopole should locate 
at every end of an open string such that ^2 V M J S/J = 
Qx 7^ 0. Condensation of these monopoles may imply 
sufficiently large fluctuations of A em and drive the system 
into the confinement phased. In fact, in the 3D compact 
case, the system is known to stay always in the confine- 
ment phas o 15 ' 16 . In the 4D compact case, there is a gauge 
transition from the confinement phase to the Coulomb 
phase as 1/e 2 increases and a Higgs transition from the 
Coulomp phase to the Higgs phase as K' increases^. 
Here let us comment on the 3D multi-component Higgs 
model in the compact version with TV Higgs fields in the 
London limits. In contrast to the above 3D model with 
N = 1, the model with N > 2 supports the Higgs phase 
due to the extra phase degrees of freedom that are free 
from coupling to the gauge field. 

The above discussion clearly shows that the SC phase 
transition do take place even in the London limit for the 
gauge Higgs model in the noncompact gauge version with 
N = 1 and in the compact version with N < 2, and 
topological excitations of SC order parameters, vortices, 
play an important role for that. 

Usually, the genuine transition temperature T c of the 
system (|2.1[) is lower than the bare critical temperature 
T° due to fluctuation effect. The radial degrees of free- 
dom of ip x may certainly contribute such renormalization 
of T c , but should not change the universality class of the 
continous phase transition we are to study because their 
fluctuations are massive. 

In the rest of the present paper, we shall study the 
FMSC state by starting with a lattice model correspond- 
ing to Eq. (|2.2p in the London limit, in which vortices are 
expected to be generated spontaneously, and therefore 
nonperturbative study is indispensable for the investiga- 
tion. 



B. GL thory in the continuum 

In the proposed GL theory^ for the FMSC materials 
in the 3D continuum space at finite T, the free energy 
density /gl for the SC state measured from the normal 
state is given as 

/gl - KJ2( D JT ■ ( D J) + a ( T - T scW ■ $ 

+A(V> ■ iff + K' ^(c^m) 2 + (T - T° m )to 2 

+«/(to 2 ) 2 + / z , 

0^ = 0^- 2ieA fll f z = -Jm ■ S, 

S = i^*x^, (2.11) 

with the spatial direction index // = 1,2,3. The three- 
component complex field ip = ("01; "02, "03)* is the spin- 
triplet SC order parameter, i.e., the Cooper-pair field (we 
omit the spatial coordinate x in the field ?p(x), etc.). tp 
is proprotional to the c?-vector in the spin space, ip oc d, 
and also the wave function of the p-wave SC state in the 
real space as a result of the spin-orbit coupling. In terms 
of ip, the magnetization (spin and angular momentum) 
S of Cooper pairs is expressed as in Eq. (|2.11[ ) 5 . The 
FM order is described by the magnetization field fh of 
electrons that do not participate in the SC state. Xg C 
and Tp M are the bare critical temperatures of SC and 
FM transitions, respectively. Because to satisfies V-m — 
0, it can be expressed in terms of the vector potential 
A as m = rotA Because the Cooper-pair field bears 
the electric charge — 2e, it couples with A minimally via 
the covariant derivative reflecting the electromagnetic 
gauge invariance. We note that this vector potential A 
is not the one for the external electro-magnetic field A%™ 
in Eq. (|2.2p . fz is the Zeeman coupling term between 
S and to. It may induce the FMSC stated, in which 
(to) ^ 0, (S) ? 0. 

The GL theory (12. lip and the related ones have been 
studied so far by means of MFT-type methods^. How- 
ever, the gauge coupling between A and makes a simple 
MF approximation assuming, e.g., a constant SC order 
parameter and ignoring topologically nontrivial fluctua- 
tions unreliable for the description of the FMSC state. 
This is also indicated as we explained for the model of 
s-w ave S C state, Eqs.([2T]) and ([2J2]). Therefore / G l in 
Eq. (|2.1ip is a kind of frustrated system of the AF and 
SC. We dare to use the word "frustrated" here. This 
is because the case with a finite magnetization to ^ 
corresponds to the case with a non-vanishing external 
magnetic field; a well known case of frustration. In fact, 
the phase of matter field there should acquires a finite 
additional phase when it is rotated along a closed loop 
and so the phase is not determined uniquely except for 
an integer magnetic flux inside the loop. In other words, 
vortices are to be generated because of the presence of 
the nonvanishing magnetization. 
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It is of course important to study the set of GL equa- 
tions derived from / GL in Eg. (12.111) . The GL equations 
should be solved self-consistently to obtain ip(x) and 
A^x). In the FM phase, ip(x) describes multi- vortex 
states and the vector potential A fJ- (x) is determined by 
the locations of vortices in addition to the other terms in 
/gl including the magnetization m(x). 

In the present paper, we shall introduce a GL theory 
defined on a cubic lattice which is a discretized version 
of Eq. (|2.1ip , and study its phase structure etc by means 
of the Monte-Carlo simulations. In this approach, all 
relevant fluctuations of tp and A^ arc taken into account. 
Some related lattice model describing a two-component 
SC was studied and interesting results were obtained 19 . 
The above two approaches, MFT of the GL solutions and 
the numerical study of the GL theory, are complementary 
and not exclusive each other. 



C. Derivation of the lattice model 

As announced in Sect. I, we introduce an effective lat- 
tice gauge model based on the GL theory in the contin- 
uum (|2.1ip by making a couple of simplifications. We 
stress that topological defects such as vortices are al- 
lowed on the lattice without introducing an additional 
short-distance cutoff for vortex cores. This point is quite 
important because it is expected that such nontrivial ex- 
citations are generated in the FMSC phase. 

As the first step of simplification, we consider the "Lon- 
don limit" of tp such that ip* -ip = const., neglecting its ra- 
dial fluctuations. As discussed in Sect. II. A for the s-wave 
SC model, this is legitimate because the phase degrees of 
freedom of i/j play an essential role for (in) stability of the 
SC state. 

Second, the FMSC materials have a FM easy-axis, 
which we choose the ^-direction (/i — 3) 20 . The Zee- 
man coupling fz prefers such that iptt x + ^2 or 
tpll oc tpi — iip2 takes large amplitude. In fact, as the 
Fermi surfaces of up and down-spin electrons have dif- 
ferent energies due to fz, the Cooper-pair amplitude of 
mixed spins, -03 = V'tii 1S small. 

Therefore, in the effective model, we ignore tps as in 
the previous studies^i 2 ^, and consider only the two- 
component complex field (tpi,ip2) in the London limit, 

(^1, 1P2, 1P3) -> (01,-02,0) 

with |Vi| 2 + IV2I 2 = ^(T s ° c -T). (2.12) 

Here we introduce a two-component complex field z that 
is the normalized Cooper-pair field; z = {zi,z 2 ) t satisfy- 
ing 

M 2 + N 2 = i, (2.13) 

and from fl2T2j) 



(^2f = J^( T ic-T)(zuZ 2 )\ (2.14) 



We note that two-component complex variables satisfy- 
ing Eq. (|2.13l) such as z(x) is called a CP 1 (complex pro- 
jective) field. In term of z(x), the first term of /gl in 
Eq. (|2.1ip is given as 

^(Tl c -T)J2(.D,zy-(D,z). (2.15) 

Now let us consider the effective lattice model on the 
3D cubic lattice. Its GL free-energy density f x at the site 
x is given up to an irrelevant constant by 

3 2 

fx = -y (Zx+n,aU X n,z xa + c.c.) - c 2 m 2 

f_L—l a—1 

-c 3 m x ■ S x + c 4 (ml) 2 -c 5 ^rh x ■ rh x+ll , 

U xfl = exp(iA xfl ). (2.16) 

The five coefficients (i = 1 ~ 5) in f|2 . 1 6[) are 
real nonnegative parameters that are expected to dis- 
tinguish various materials in various environments. z x = 
(zxijZ^Y Q2 a= \z xa \ 2 — f) is the CP 1 variable at the 
site x and plays the role of SC order-parameter field. 
U xfl is the exponentiated vector potential^! put on the 
link (x,x + /i). rh x = (m x i, m X 2, w X 3) t is the magnetic 
field made out of A xfl as 

3 

m ¥ = c^aV^Aea, V^Aea = A x+Vt x-A x \. (2.17) 

i/,A=l 

rh x serves as the FM order-parameter field. S x = 
(0, 0, S x sY is an Ising-type spin vector of Cooper pairs 
made out of z xa as 

5x3 = i[z xl z x2 - z* 2 z xl ) oc \ip n \ 2 - |0u| 2 - (2-18) 

As for the case of (|2.7p . f x is invariant under a U(l) gauge 
transformation, 

Z X a ^ Z xa CXp(z ^)Z xa , 

U xfl -> U' xfi = cxp(i\ x+tl )U xfi exp(-£A x ). (2.19) 

We note that both rh x and S x are gauge invariant. 

The meaning of e£Lch term in f x is as follows. The ex- 
terna describes a hopping of Cooper pairs. From (12.151) . 
it is obvious that 

Cl ~ A"aA _1 (Tg C - T)a, (2.20) 

where a is the lattice spacing. At sufficiently large c\, 
the ci-term may stabilize the phase of z xa , and then a 
coherent condensation of the phase degrees of freedom of 
z x induces the superconductivity. The C2 and C4-terms 
are the quartic GL potential of rh x and favor a finite 
amount of local magnetization (rh x ) =/= (note that we 
take C2 > 0). Again we note these terms controls in- 
trinsic magnetrization and different from Fa of (|2.2p 
for the fluctuationg but external magnetic field. The 
cs-term enhances uniform configurations of m x , i.e., a 
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FM long-range order signaled by a finite magnetization, 
limi x _, r /|_ > . 00 (TO a : ' ™>x') 7^ 0- The C3-term is the Zeeman 
coupling, which favors collinear configurations of rh x and 
S x , namely the coexistence of ferromagnetism and super- 
conductivity. 

The partition function Z at T is given by the integral 
over a set of two fundamental fields z xa and A x/Jl as 

Z= [ [dz][dA] exp(-/?n (3 = T-\ F = £/ B , 
[dz] = \\d 2 : z xl d 2 z x2 S(\z xl \ 2 + \z x2 \ 2 - 1), 

X 

[dA] = Y[ dA xt ,, A XI1 g (-co, oo). (2.21) 

The coefficients Cj in may have nontrivial In- 
dependence as Eas. (|2.1ip and (12.20)) suggest. However, in 
the present study we consider the response of the system 
by varying the "temperature" T = 1//3 defined by /?, an 
overall prefactor in Eq. (|2.21[ ). while keeping Cj fixed. This 
method corresponds to well-known studies such as the 
FM transition by means of the 0(3) nonlinear-cr model 22 
and the lattice gauge- Higgs models discussed in Sect. II. A, 
and sufficient to determine the critical temperature. See 
later discussion leading to Eq. (l3.6[) . 

III. NUMERICAL STUDIES 



A. FM ands SC phase transition and Meissner 
effect 



For explicit procedures of our Monte-Carlo simula- 
tions, we first prepare a 3D lattice of the size of (2 + 
L + 2) 2 x L, namely the lattice coordinates runnning as 
x%, X2 = 0, L+3, x 3 = 0, L— 1. The extra width 2+2 
in the \x — 1,2 directions is introduced as a buffer zone 
in which the suppercurrrent is damped. The calculations 
of physical quantities are done in the central region R of 
the size L 3 to suppress the effects of the boundary. For 
the boundary condition we choose the periodic boundary 
condition in the fi = 3 direction. For the fx = 1, 2 direc- 
tions, we first note that the supercurrent density on 
the lattice is expressed as 



y -"X J • 



(3.1) 



Then we impose 



Zx+fi - U xfl z x = for x = (0, x 2 , x 3 ), fj, = I, 

x = (L + 2,x 2 ,x 3 ),fi = 1, 

x = (xi, 0, £3), /i = 2, 

x = (xi,L + 2,x 3 ),fi = 2, (3.2) 

whereas we impose the free boundary condition on A Xil . 

The above condition (I3.2[) is gauge invariant and as- 
sures us that = (fi — 1, 2) on the boundary surfaces 
in the fx— 3 plane; the supercurrent do not leak out of the 



SC material put in the region R. We note this boundary 
condition enhances the third-component of the FM order 
so that (rn x ) = (0, 0, m^f . This fact is traced back to our 
choice 1P3 — in Eq. (|2.12p and reflects the experimental 
fact that the real FMSC materials exhibit the Ising-type 
anisotoropy of the FM magnetization. Note also that 
the conventional periodic boundary condition for A x/X in 
all the three directions implies J2 X = an d the FM 
order cannot exist. 

We use standard Metropolis algorithm for the lattice 
size up to L = 20. The typical sweeps for measurement is 
(30000 - 50000) x 10 and the acceptance ratio is 40% ~ 
50%. 

Explicitly, we calculate the internal energy U, the spe- 
cific heat C of the central region R, the FM magnetiza- 
tion m^, 

U = ^(F'}, C = -L<(f - (F'}) 2 >, F' = £ f x , 



(3.3) 



xG-R 



and the normalized correlation functions, 



G m (x - x ) = (m x ■ rn xo )/(m xo ■ m xo ), 

G s (x - x ) = (S x3 S XOy3 )/(S XOy3 S XOy3 ), (3.4) 

where xq is chosen on the boundary of R such as (3, 2 + 
L/2,z). 

We first show that the model ([2T2T1) exhibits a FM 
phase transition as T is lowered. To this end, we put 
ci — c 3 — and (02,04,05) = (0.5,4.0,1.0), and increase 
j3 in the Boltzmann factor of Eq. (|2.21j ). In FigO]we show 
C and m M . It is obvious that a second-order phase transi- 
tion to the FM state takes place at /3fm = 1/7fm — 2.0. 
We verified that other cases with various values of 02,4,5 
exhibit similar FM phase transitions. 

Next, let us study the SC phase transition. Here it 
is useful to consider the case of all Ci — except for 
ci. Then the model is related to the CP 1 +U(1) lattice 
gauge theory 2 ^ which has the energy of the form f r = 
~(c\/2)^(zUz + c.c.) + Fa where Fa is the compact 
version (|2.5p . In fact, they agree by setting 1/e 2 = 0. 
The phase structure of this model is studied in Ref. 23 and 




3.5 4 



FIG. 1. (a) Specific heat for (c 2 ,c 4 ,e 5 ) = (0.5,4.0,1.0) and 
ci = C3 = 0. At j3 ~ 2.0, C exhibits a sharp peak indicating 
a second-order FM phase transition, (b) Each component of 
magnetization m M vs j3. For T < Tfm, mz develops, whereas 
mi and m.2 are zero within the errors as expected. 
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FIG. 2. (a) Specific heat vs /3 for (ci,c 3 ) = (0.2,0.2) and 
(c2,C4,c 5 ) = (0.5,4.0,1.0) (L = 20). There are a large peak 
at P ~ 2.1 and a small one at /3 ~ 4.5. (b,c) Specific heat d 
of each term of F in (|2. 2 1 f) vs /3. The small and broad peak 
at /? ~ 4.5 in C is related to fluctuations of ci-term. 



it is found that the phase transition from the confinnicnt 
phase to the Higgs phase takes place at c\ ~ 2.85 for 
1/e 2 = 0. Thus the SC state exists at sufficiently large 
c x . 

For simulation, we put (01,03) = (0.2,0.2) and 
(c 2 ,c 4 ,c 5 ) = (0.5,4.0,1.0). In FiglSi, we show C vs /3. 
There are a large and sharp peak at j3 ~ 2.1 and a small 
and broad one at /3 ~ 4.5. In order to understand phys- 
ical meaning of the second peak, it is useful to measure 
"specific heat" of each term fa in the free energy (|2.16p 
defined by C t = ((F- - (Fl)) 2 )/L 3 . FigHt shows that 
the specific heat of the ci-term has a relatively large and 
broad peak at j3 ~ 4.5. Then we conclude that the SC 
phase transition takes place at /?sc — 4.5. 

To verify this conclusion, we show G m (r) and Gs(r) 
in FigO At j3 = 2.5, G m (r) exhibits a finite amount 
of the FM order, whereas Gs(r) decreases very rapidly 
to vanish. This means that, as T is decreased, the FM 
transition takes place first and then the SC transition 
does. Therefore, for /3 > /3 S c - 4.5, the FM and SC 
orders coexist. 

It is interesting to clarify the relation between the bare 
transition temperature Tg C in Eq. (|2.1ip and the genuine 
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FIG. 3. Correlation functions G m (r) and Gs(r) at various 
T's for L = 20. Cj's are the same as in FigO 



transition temperature T$c- From Eq. (|2.2ip . any physi- 
cal quantity is a function of /3c,. In the numerical simu- 
lations, we fix the values of Ci and vary j3 as explained. 
Then the result /3sc — 4.5 means 



/3ci|t=t sc = 4.5 x 0.2. 



(3.5) 



By using Eq. (l2.20p . this gives the following relation; 



1 Ka(T° c - T sc )a 



Tsc 
Tsc — 



A 

0.90 A 

Ka a 



= 0.90, 



sc- 



(3.6) 



Eq. p.6p shows that the transition temperature is lowered 
by the fluctuations of the phase degrees of freedom of 
Cooper pairs. We expect that a relevant contribution 
to lowering the SC transition temperature comes from 
vortices that are generated spontaneously in the FMSC 
as we show in Sec. III. B. 

After having confirmed that the genuine critical tem- 
perature can be calculated by the critical value of ft with 
fixed Ci, we use the word temperature in the rest of the 
paper just as the one defined by T = 1//3 while Ci are 
T-independent parameters. 

Meissner effect is one of the most important phenom- 
ena for a SC order. To study it, we follow the following 
steps2£; (i) introduce a vector potential A™ for an ex- 
ternal magnetic field, (ii) couple it to Cooper pairs by 
replacing U xfl — >■ U x)J ,exp(iA'£l l ) in the c\ term of f x and 
add its magnetic term /£ x = +C2'(m® x ) 2 (c' 2 > 0) to 
f x with m® x defined in the same way as (|2.17p by using 
A° x , (iii) let A™ fluctuate together with z xa and A xfx 
and measure an effective mass Mq of A° x via the decay 
of correlation functions of m™. The result of m™ prop- 
agating in the 1-2 plane is shown in FigfJ] It is obvious 
that the mass Mq starts to develop at the SC phase tran- 
sition point, and we conclude that Meissner effect takes 
place in the SC state. 
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FIG. 4. Gauge-boson mass Mg of the external magnetic field 
propagating in x-y plane vs f3 for the same a as in Fig[2] 
and C2' = 3.0 (L — 16). At the SC phase transition point 
/3sc — 4.5 (indicated by an arrow) determined by the peak of 
C, Mq starts to increase from small values. 
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FIG. 5. The results for a constant magnetic field fn x = 
(0,0,771^3)' at ci = C3 = 0.2 (L = f6). (a,b) Specific heat 
C and (c,d) the real part of SC correlation function G z (r) in 
the 1-2 plane. (a,c) m X 3 — n/4 and (b,d) m X 3 = 7r. 



where mod(a;) = mod(:r, 2ir), In short, V^r describes 
vortices of electron pairs with the amplitude V>-ff(m = 
ipi ± oc zi ± iz2. 

In FigJS]we present snapshots of V x at c\ = c 3 = 0.2 
for fixed values of gauge potential A xfl corresponding to 
a constant magnetization, rh — (0,0, ?)*. It shows that 
(i) both of the fluctuations around zero, ((V^l), decrease 
as P increases, and (ii) V x has larger fluctuations than 
V~ at high T, whereas smaller ones at low T. These 
behaviors are consistent with the Zeeman C3-term in the 
energy f x of (|2.16[) . which distinguishes the z x order and 
the z~ order, and the fact that m directs to the third- 
direction in the present case. In Figj7l we also show the 
vortex snapshots at ci = c 3 = 0.2 and fa — (0,0,77)*. 
Compared with the case m — (0, 0, vortices here are 
located rather systematically as we expected from the 
result of correlation function G z {r). 



B. SC transition and vortices in a constant 
magnetic field 

Because the observed SC state in Figs|2] and [3] coexists 
with the FM order, it is expected that vortices of the SC 
order parameter are induced there spontaneouly 2 ^. To 
verify this expectation, we set the vector potential A xfl 
to a position-dependent but nonfluctuating value that 
corresponds to a uniform magnetic field in the third- 
direction, and study the behavior of z xa itself. In this 
case, the free energy f x loses the local gauge symmetry 
(|2.19j) , and therefore the correlation function of z xa , 

G z (x - X ) = (z x ■ Zxg), ( 3 - 7 ) 

has nonvanishing values in the SC phase. 

In FigJSl we show C and G z (x) for two cases of fixed 
rh x . For the case with m x — (0,0,^)', C has a shape 
similar to C\ in FigJ5J;, and indicates a SC phase tran- 
sition at p ~ 4.8. G z (r) exhibits fluctuating behavior 
even for low T's, j3 > 4.9. This suggests that vortices are 
spontaneouly generated in the SC state violating spa- 
tial unifomity, and their locations fluctuate. In the other 
case of rh x = (0, 0, it) 1 , C has a sharper peak at ft ~ 4.5, 
and G z (r) exhibits clear periodically oscillating behavior 
with the period 4 x (lattice spacing). This implies that, in 
this case, locations of vortices are rather stable compared 
with the case of m x ^ = 

In order to verify the above expectation, we calculate 
vortex density directly. In the present model, one may 
define the following two kinds of gauge-invariant vortex 
densities V x and V~ in the 1-2 plane; 



z^ = z xl ± iz x2 = 


exp(it 


^) 

'x 11 


Vt = ^[mod(,9± 


-1 " 




A xl ) 


+mod(0± +1 _ 


h2 


- 


~ A x 


-mod(0j +1 _ 






-A x 


-mod(0^ +2 


- < 







C. Order of FM and SC phase transitions and 
phase diagram 

Let us next examine the possibility that the order of 
the FM and SC phase transitions are interchanged as the 
value of C5 is decreased. As the C5-term in f x tends to 
align m x , Tfm decreases as C5 is decreased. Most of the 
FMSC materials loses the FM order as the applied P 
is increased, and then it is phenomenologically expected 
that C5 is a decreasing function of P. We consider several 




FIG. 6. Snapshots of vortex densities V x at ci = C3 = 0.2 
for a fixed m = (0,0,tt/4) (L = 16). Black dots; V x ± = 0.875, 
Dark gray dots; -1.125, Light gray dots; -0.125. (A)V+ at /3 = 
3.0, (B)KT at P = 3.0, (C)V+ at p = 7.0, (D)V7 at P = 7.0. 
The average magnitude {\V X ±\) is (A) 0.387, (B) 0.380, (C) 
0.331 and (D) 0.335. The points Vt = -0.125 = -m 3 /(2n) 
reflect m itself. 
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FIG. 7. Snapshots of vortex densities at C\ = C3 = 0.2 for 
a fixed m = (0,0, tt) (L = 16). (A)V+ at (3 = 3.0, (B)V~ at 
P = 3.0, (C)V+ at /3 = 7.0, (D)V~ at /3 = 7.0. The average 
magnitude flV^I) is (A) 0.528, (B) 0.537, (C) 0.560 and (D) 
0.508. The points Vt = -0.5 = -m 3 /(27r) reflect m itself. 



cases with C5 = 1.5, 1.0,0.7,0.5,0.4 and 0.3, while other 
Cj are fixed to the same values as those in Figj2] where 
c 5 = 1.0. 

We find that the order of two phase transitions actually 
interchanges at C5 ~ 0.5. In FigEl we show the specific 
heat G,Cx,C^ for C5 = 0.4. C has the two peaks at 
fisc ~ 5.0 < /?fm ~ 8.3. C\ is sharper than in the case 
of c 5 = 1.0 in FigH 

In Fig [9] we present G m {r) and Gs{r) with r in the 
1-2 plane for C5 = 0.4, which exhibit very peculiar be- 
havior; In the FM and SC coexisting phase (T < Tfm), 
they have nonvanishing values only near the surface of 
the lattice in contrast with FigJ3] This behavior sur- 
vives in larger systems. For example, we define the 
thickness AL of the coexisting region such that the 
ordered region in the 1-2 plane occupies the interval 
2 + AL + (disordered region) + AL + 2 in the lattice length 
2 + L + 2 in the fx = 1, 2 directions. We obtain AL ~ 3 
for L = 12(FigiJ and AL ~ 4 for L = 16, so about 
the outer half region in the linear dimension is occupied 
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FIG. 8. Specific heats for c 5 = 0.4 (L = 12). (a) Total C 
and the specific heat Ci of the Ci-term, (b) C5 of the cs-term. 
Order of two phase transitions is interchanged (Psc < /3fm)- 
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FIG. 9. Correlation functions G m (r) and Gs{r) at C5 = 0.4 
(L = 12). The other a are the same as those in Fig[2] They 
exhibit orders near the surface of the lattice. 

by the ordered state. This implies that the FMSC co- 
existing phase appear in the region including the surface 
of the material, and not in the central region inside the 
system. We note that this "surface" region is not two- 
dimensional but 3D, because the SC-FM transition is a 
genuine second-order one, which is not allowed in a two- 
dimensional system^. This phenomenon is a prediction 
of the present model. 

It is intriguing to draw a phase diagram in the P- 
T plane assuming certain phenomenological relation be- 
tween C5 and P. In the experiments, the critical tem- 
perature Tfm is a decreasing function of P. This means 
that the parameters C2, C4 and C5 in Eg. (12.161) vary as 
functions of P. Changes of C2 and/or C4 influence the 
magnitude of the magnetization vector rh x and result in 
a change of C5 after a renormalization of m x . Then for 
example, one may "phenomenologically" assume 

C5 = c 5 \l —J , (3.9) 

where P c is the critical pressure at which the FM order 
disappears even at T — > (i.e., at c 5 = 0), c\ is the 
value of C5 at which P = 0, and the power 7 is a fitting 
parameter. In Fig llOl we show the phase diagram drawn 
with certain choice of these parameters. This phase dia- 
gram has a similar structure to the experimental result 
of UCoGe. 

In the same way as decreasing C5, the case of increas- 
ing C3 has been studied with several choices of C5 and 
(ci,C2,C4) = (0.2,0.5,4.0). Both Tfm and Tsc increase 
as C3 increases. Furthermore, for sufficiently large val- 
ues such as C3 = 1.5, Tfm > Tsc even for C5 = 0.5 as 




0.5 1 c 5 1.5 2 0.2 0.4 0.6 0.8 1 P(C 5 ) 



FIG. 10. Phase diagrams in (a) C5-T and (b) P-T planes 
with P c = 1.0, C5 = 2.2 and 7 = 3.0. ci — C4 are same as in 
Fig|2] 
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expected. This indicates that the present model with 
larger C3 may provide a phase diagram similar to that of 
UGe 2 and URhGe. 

IV. CONCLUSION 

In summary, we have proposed a GL model defined on 
the 3D lattice for the FMSC state, and shown that it 
explains some experimental observations such as the 
phase diagram and the homogeneous and inhomogeneous 
FMSC states. This model naturally includes effects of 
topological excitations, vortices, that play an important 
role for the SC phase transition in the FM state. Al- 
though the obtained global phase structure is similar to 
that of MFT, the appearance of inhomogeneous configu- 
rations such as the FMSC state and vortex configurations 
are certainly beyond the scope of MFT. In the present 
analysis, we consider the "London limit", in which the 



radial fluctuations of the two-gap SC order parameters 
are ignored. As we explained in Sec. II. A., these fluc- 
tuations may change the order of SC phase transitions 
and may play an important role in SC transitions that 
are induced by an external magnetic-field. This problem 
is under study and results will be reported in a future 
publication. 
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